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Abstract. We provide a rigorous Bayesian formulation of the EIT problem 
in an infinite dimensional setting, leading to well-posedness in the Hellinger 
metric with respect to the data. We focus particularly on the reconstruction 
of binary fields where the interface between different media is the primary 
unknown. We consider three different prior models - log-Gaussian, star-shaped 
and level set. Numerical simulations based on the implementation of MCMC 
are performed, illustrating the advantages and disadvantages of each type of 
prior in the reconstruction, in the case where the true conductivity is a binary 
field, and exhibiting the properties of the resulting posterior distribution. 


1. Introduction 

1.1. Background. Electrical Impedance Tomography (EIT) is an imaging tech¬ 
nique in which the conductivity of a body is inferred from electrode measurements 
on its surface. Examples include medical imaging, where the methodology is used 
to non-invasively detect abnormal tissue within a patient, and subsurface imaging 
where material properties of the subsurface are determined from surface (or occa¬ 
sional interior) measurements of the electrical response; the methodology is often 
referred to as electrical resistance tomography - ERT - in this context and discussed 
in more detail below. The concept of EIT appears as early as the late 1970’s [15] 
and ERT the 1930’s [28]. 

A very influential mathematical formulation of the inverse problem associated 
with EIT dates back to 1980, due to Calderon. He formulated an abstract version 
of the problem, in which the objective is recovery of the coefficient of a diver¬ 
gence form elliptic PDE from knowledge of its Neumann-to-Dirichlet or Dirichlet- 
to-Neumann operator. Specifically, in the Dirichlet-to-Neumann case, if D C 
and g G is given, let u G H^{D) solve 

J— V • {a\/u) = 0 in D 

1 u = g on dD. 
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The problem of interest is then to ask does the mapping A^- : (dD) 

given by 


du 

determine the coefficient a in D1 Physically, g corresponds to boundary volt¬ 
age measurements, and Acr{g) corresponds to the current density on the boundary. 
Much study has been carried out on this problem - some significant results, in the 
case where all conductivities are in C‘^{D) and d > 3, concern uniqueness [42], 
reconstruction [35], stability [2] and partial data [22]. Details of these results are 
summarised in [38]. 

In 1996, Nachman proved global uniqueness and provided a reconstruction pro¬ 
cedure for the case d = 2, involving the use of a scattering transform and solving 
a D-bar problem [36]. The D-bar equation involved is a differential equation of 
the form = /, where d denotes the conjugate of the complex derivative and / 
depends on the scattering transform. A regularised D-bar approach, involving the 
truncation of the scattering transform, was provided in [23, 24], enabling the recov¬ 
ery of features of discontinuous permeabilities. The regularised D-bar approach is 
also used in [25], for the case when the data is not of infinite precision. Other work 
in the area includes joint inference of the shape of the domain and conductivity [26]. 

For applications, a more physically appropriate model for FIT was provided in 
1992 in [40]. This model, referred to as the complete electrode model (GEM), re¬ 
places complete boundary potential measurements with measurements of constant 
potential along electrodes on the boundary, subject to contact impedances. The au¬ 
thors show that predictions from this model agree with experimental measurements 
up to the measurement precision of 0.1%. For this model they also prove existence 
and uniqueness of the associated electric potential. It is this model that we shall 
consider in this paper, and it is outlined in section 2. 

When using the GEM, there is a limitation on the number of measurements that 
can be taken to provide additional information. This is due to the fact that there 
are only a finite number of electrodes through which current can be injected and 
voltages read, combined with the linear relationship between current and voltage. 
The data is therefore finite dimensional in the inverse problem, as distinct from the 
Galderon problem where knowledge of an infinite dimensional operator is assumed. 
As a consequence, reconstruction using the GEM often makes use of Tikhonov 
regularisation; such a regularisation can also be used to account for noise on the 
data. The paper [13] analyses numerical convergence when an or TV penalty 
term is used, with a finite element discretisation of the problem. We will adopt a 
Bayesian approach to regularisation, and this is discussed below. 

A closely related problem to FIT is Electrical Resistivity Tomography (ERT), 
which concerns subsurface inference from surface potential measurements, see for 
example [28] which discussed the problem as early as 1933. Physically the main 
difference between FIT and ERT is that alternating current is typically used for the 
former, and direct current for the latter. Additionally, due to the scale of ERT, it 
is a reasonable approximation to model the electrodes as points, rather than using 
the GEM. Another difference between the two is that the relative contrast between 
the conductivities of different media are typically higher in subsurface applications 
than medical applications, which permits the approximation of the problem by 
a network of resistors in some cases [37]. Nonetheless, the Bayesian theory and 
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MCMC methodology introduced here will be useful for the ERT problem as well as 
the EIT problem. 

Statistical, in particular Bayesian, approaches to EIT inversion have previously 
been studied, for example in [18, 19, 31]. In [18], the authors prove certain regu¬ 
larity of the forward map associated with the CEM, formulate the Bayesian inverse 
problem in terms of the discretised model, and investigate the effect of different 
priors on reconstruction and behaviour of the posterior. The paper [31] focuses on 
Whittle-Matern priors, using EIT and ERT as examples for numerical simulation. 
The paper [19] presents a regularised version of the inverse problem, which admits 
a Bayesian interpretation. 

The Bayesian approach to inverse problems, especially in infinite dimensions, 
is a relatively new technique. Two approaches are typically taken: discretise first 
and use finite dimensional Bayesian techniques, or apply the Bayesian methodology 
directly on function space before discretising. The former approach is outlined in 
[20]. The latter approach for linear problems has been studied in [12, 32, 33, 34]. 
More recently, this approach has been applied to nonlinear problems [10, 29, 30, 41]. 
It is this approach that we will be taking in this paper. 

1.2. Our Contribution. The main contributions in our paper are as follows: 

(i) This is the first paper to give a rigorous Bayesian formulation of EIT in infinite 
dimensions. 

(ii) This setting leads to proof that the posterior measure is Lipschitz in the data, 
with respect to the Bellinger metric, for all three priors studied; further stabil¬ 
ity properties of the posterior with respect to perturbations, such as numerical 
approximation of the forward model, may be proved similarly. 

(hi) We employ a variety of prior models, based on the assumption that the 
underlying conductivity we wish to recover is binary. We initially look at 
log-Gaussian priors, before focusing on priors which enforce the binary field 
property. These binary field priors include both single star-shaped inclusions 
parametrised by their centre and by a radial function [8], and arbitrary geo¬ 
metric interfaces between the two conductivity values defined via level set 
functions [17]. 

(iv) Numerical results using state of the art MCMC demonstrate the importance 
of the prior choice for accurate reconstruction in the severely underdetermined 
inverse problems arising in EIT. 

1.3. Organisation of the Paper. In section 2 we describe the forward map asso¬ 
ciated with the EIT problem, and prove relevant regularity properties. In section 
3 we formulate the inverse problem rigorously and describe our three prior models. 
We then prove existence and well-posedness of the posterior distribution for each of 
these choices of prior. In section 4 we present results of numerical MCMC simula¬ 
tions to investigate the effect of the choice of prior on the recovery of certain binary 
conductivity fields. We conclude in section 5. 

2. The Eorward Model 

In subsection 2.1 we describe the complete electrode model for EIT as given in 
[40]. In subsection 2.2 we give the weak formulation of this model, stating assump¬ 
tions required for the quoted existence and uniqueness result. Then in subsection 
2.3 we define the forward map in terms of this model, and prove that this map is 
continuous with respect to both uniform convergence and convergence in measure. 
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2.1. Problem Statement. Let I) C (i < 3, be a bounded open set with 
smooth boundary representing a body, with conductivity a : D A number L of 

electrodes are attached to the surface of the body. We treat these as subsets 
of the boundary dD, and assume that they have contact impedances G 

A current stimulation pattern G is applied to the electrodes. Then the 

electric potential v within the body and boundary voltages (V/)^i ^ 
are modelled by the following partial differential equation. 


( 1 ) 


—V • (a(x)Vv(x)) =0 X e D 


L 


a^dS = Ii 

on 

/ N / N 

v{x) + Zia{x)^{x) = Vi 


1 = 1,...,L 
X €dD\ e; 


X G e/, / = 1,..., 1/ 


This model was first proposed in [40]; a derivation can be found therein. Note 
that the inputs to this forward model are the conductivity a, input current 
and contact impedances The solution comprises the function v : D R 

and the vector V G R^ of voltages. Also note that solutions to this equation are 
only defined up to addition of a constant: if (v^V) solves the equation, then so does 
(n + c, y + c) for any c G M. This is because it is necessary to choose a reference 
ground voltage. 



Figure 1. An example domain D with electrodes attached 

to its boundary. 
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2.2. Weak Formulation. We first define the space in which the solution to equa¬ 
tion (1) will live. Using the notation of [40] we set 

\\iv,V)\\l=\\v\\l.^j,) + \\V\\% 

= \\v\\UD) + \\^v\\h^D) + \\V\\%. 

Since solutions are only defined up to addition of a constant, we define the quotient 
space (HI, || • ||^) by 


HI = H/R, 

||(t;,y)||e = inf ||(t,-c,y-c)||H. 

cGM 

We will often use the notation v' = (v^V) for brevity. It is more convenient to equip 
H with an equivalent norm, as stated in the following lemma from [40]: 

Lemma 2.1. Define ||'r’'||* by 

1=1 

Then || • ||>^ and || • ||j^ are equivalent. 

We can now state the weak formulation of the problem as derived in [40]. For 
this let w' = (re, W). 

Proposition 2.2. Let 5 : El x El ^ M and r : El ^ M be defined by 

B{v',w';a)= [ crVi; ■ Vw da; + V 1 / {v - Vi){w - Wi) dS, 

Jd Zl Jei 

L 

r{w') = J2liWi. 

1 = 1 

Then if v' is a strong solution to the problem (1), it satisfies 
(2) B{v'a) = r{w') forallw'eM. 


We will use the weak formulation (2) to define our forward map for the complete 
electrode model (1). In order to guarantee a solution to this problem, we make the 
following assumptions. 

Assumptions 2.3. The conductivity a : I) ^ M, the contact impedances (zi)fLi G 
and the current stimulation pattern ^ satisfy 

(i) a G L^{D;R), essinf a{x) = cr_ > 0; 

x^D 

(ii) 0 < Z- < zi < < oc, / = 

L 

(iii) J2li=0. 

1 = 1 

Under these assumptions, existence of a unique solution to (2) is proved in [40] 
and stated here for convenience: 
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Proposition 2.4. Let Assumptions 2.3 hold, then there is a unique [(p, G El 
solving (2). We may, without loss of generality, ehoose the element {v, V) G [(p, V)] 
of the equivalenee elass of solutions to he that whieh satisfies 

(3) E^*=0- 

1=1 

Remark 2.5. Assumptions 2.3 (i) and (ii) are to ensure coercivity and boundedness 
of •;cr). Assumption 2.3 (iii) is necessary for continuity of r(-), and physically 
may be thought of as a conservation of charge condition. Choosing a solution from 
the equivalence class corresponds to choosing a reference ground voltage. 

2.3. Continuity of the Forward Map. In what follows we will restrict to the set 
of admissible conductivities, which is defined as follows. 

Definition 2.6. A conductivity field a : D ^ M is said to be admissible if 

(i) there exists A" G N, open disjoint subsets of D for which D = 

(ii) cr|^ G C{Dj)] and 

(iii) there exist a~ G (0, oc) such that a~ < a(x) < for all x G D. 

The set of all such conductivities will be denoted A{D). 

Note that any a G A{D) will satisfy Assumptions 2.3(i). Assume that the current 
stimulation pattern G and contact impedances G are known 

and satisfy Assumptions 2.3(ii)-(iii). Then we may define the solution map Ai : 
A{D) ^ El to be the unique solution to (2) satisfying (3). The above existence and 
uniqueness result tells us that this map is well-defined. 

In [18] it is shown that M. : A{D) ^ El is Frechet differentiable when we equip 
A{D) with the supremum norm. Though this is a strong result, this choice of norm 
is not appropriate for all of the conductivities that we will be considering. We hence 
establish the following continuity result. 

Proposition 2.7. Fix a eurrent stimulation pattern G and eontaet 

impedanees (^z)/li G satisfying Assumptions 2.3. Define the solution map 
A4 : A(D) ^ El a5 above. Let a G A{D) and let C A{D) be sueh that 

either 

(i) (7^ eonverges to a uniformly; or 

(ii) (7^ eonverges to a in (Lebesgue) measure, and there exist G (0, oo) sueh 

that for all e > 0 and x e D, a~ < (7s{x) < cr+. 

Then \\M{(7A) - Af(cr)||* ^ 0. 

Proof. Define the maps 5 : El x El x A{D) M and r : El ^ M as in Lemma 2.2, 
but on El rather than El. Then denoting u'^ = {v^,V^) = Ad{(7A) and v' = = 

Af (cr), we have for all w' G El, 


B{v(, w'; cTg) = r(w'), B{v', w'; a) = r{w'). 
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It follows that 


0 = cr^) - B{v', w'; + B{v', w'; - B{v', w'; a) 

= [ a,Viv,-v)-Vwdx + y2 - [ iiv,-v)-{Vt-Vi))iw-Wi)dS 

JD l=\ 

+ / (cTg — cr)V'i; • Vr; dx. 

Jd 

Letting w' = {v^ — — V)^ we see that 

f cT,|V(?;,-ei"dx +Y- [ i{v,-v)-{Vt-Vi)fdS 

J D l=\ 

< / \(7s — (j\\S/v'S/{vs — v)\dx. 

Jd 

In both cases (i) and (ii), we have that is bounded uniformly below by 

a positive constant. Hence for small enough 5 , the left hand side above may be 
bounded below by C\\v'^ — v'\\^. We then have by Cauchy-Schwarz 


Wvi-v'Wl 


( 4 ) 

( 5 ) 


<C / \(Ts — (t\\S/v— v)\dx 

Jd 

“^(/ • l|V(ve -i;)||l2 

<C^y \(Je - (j\^\Vv\^ dx^ -Wv'e-v'W* 

< C||(Te -C^||oo||Vv||i2||vY^''ll*• 


If (Je ^ cr uniformly, we deduce from (5) that \\v'^ ~ '^^11* result follows. 

If IcTg — cr| -> 0 in measure, then since \D\ < 00 , it follows that the integrand in (4) 
tends to zero in measure, see for example Corollary 2.2.6 in [6]. Since is assumed 
to be uniformly bounded, the integrand is dominated by a scalar multiple of the 
integrable function jVt’p. We claim that this implies that the integrand tends to 
zero in L^. Suppose not, and denote the integrand /g. Then there exists (5 > 0 
and a subsequence {fsi)i>i such that ^ ^ for all i. This subsequence still 

converges to zero in measure, and so admits a further subsequence that converges 
to zero almost surely. An application of the dominated convergence theorem leads 
to a contradiction, hence we deduce that tends to zero in and the result 
follows. □ 


Denote the projection H : H ^ M^, {v,V) i-A V. The following lemma shows 
that the above result still holds if we replace Ai by H o Al. 

Corollary 2.8. Let the assumptions of Proposition 2.7 hold. Then 

|n o —Ho Ai(a )\£2 0. 

Proof. We show that there exists C > 0 such that for all {v,V) G HI with — 

0, IK'^’, > C\V\£ 2 . By the equivalence of || • and || ' IIei’ Lemma 2.1, we have 
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The infimum on the right-hand side is attained at 

1=1 

Then by Proposition 2.7, we have 

0< \IloM{(j,)-IloM{(j)y < 


□ 


3. The Inverse Problem 

We are interested in the inverse problem of determining the conductivity field 
from measurements of the voltages on the boundary, for a variety of input 

currents on the boundary. To this end we introduce the following version 

of Ohm’s law. Observe that the mapping / i-A p', taking the current stimulation 
pattern to the solution of (2), is linear. Then given a conductivity field a G A{D)^ 
there exists a resistivity matrix R{(j) G such that the boundary voltage 

measurements V(a) arising from the solution of the forward model are related to / 
via 

V(a) = R(a)I 

By applying several different current stimulation patterns we should be able to infer 
more about the conductivity a. Note however that since the mapping I V is lin¬ 
ear, only linearly independent stimulation patterns will provide more information^. 
Since we have the conservation of charge condition on /, there are at most L — 1 
linearly independent patterns we can use. 

Assume that J linearly independent current patterns G j = 1,..., J, 
J < L — 1 are applied, and noisy measurements of = R{a)I are made: 

Uj = v^^'> + ~ A^(0, To) iid. 

We have 

Vj = + »7i 

where ^j(cr) = Concatenating these observations, we write 

y = G{cr) + ri, ri'^N{0,T) 

where P = diag(ro,... ,ro). The inverse problem is then to recover the conduc¬ 
tivity field a from the data y. This problem is highly ill-posed: the data is finite 
dimensional, yet we wish to recover a function which, typically, lies in an infinite 
dimensional space. We take a Bayesian approach by placing a prior distribution on 
cr. The choice of prior may have significant effect on the resulting posterior distri¬ 
bution, and different choices of prior may be more appropriate depending upon the 
prior knowledge of the particular experimental set-up under consideration. 

In subsection 3.1 we outline three different families of prior models, and show 
the appropriate regularity of the forward maps arising from them. In subsection 
3.2 we describe the likelihood and posterior distribution formally, before rigorously 
proving that the posterior distribution exists and is Lipschitz with respect to the 
data in the Hellinger metric. 

there is noise on the measurements, additional linearly dependent observations can be made 
to effectively reduce the noise level on the original measurements. We can assume that this has 
been done and scale the noise appropriately. 
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3.1. Choices of Prior. In this section we consider three priors, labelled by i = 

1, 2, 3, defined by functions Fi : Xi ^ which map draws from prior measures 

on the Banach spaces Xi to the space of conductivities A{D). Our prior conductiv¬ 
ity distributions will then be the pushfoward of the prior measures by these maps 
Fi. We describe these maps, and establish continuity properties of them needed 
for the study of the posterior later. In what follows, the space Cq{D) of contin¬ 
uous functions on D will be equipped with the supremum norm II ■ Iloo and the 
corresponding Borel cr-algebra. 

3.1.1. Log-Gaussian prior. We first consider the simple case that the coefficient is 

given by the exponential of a continuous function. Let Fi : C^{D) -4(D) be 

defined by Fi{u) = exp(R). Then it is easily seen that Fi does indeed map into 
A{D). Furthermore, since D is bounded, if r G C^{D) and {Ue)e>0 C C°{D) is a 
sequence such that \\us — r||oo ^ 0, then \\Fi{us) — Di(r)||oo ^ 0- 

In this case, we will take our prior measure /io on r to be a non-degenerate 
Gaussian measure N{mo,Co) on C^{D). Note that the push forward of a Gaussian 
measure by Fi is a log-Gaussian measure. 

Example 3.1. Gonsider the case D = D(0,1) C Suppose that u is drawn 
from a Gaussian measure /io = A/'(0,C). Typical samples from Ff(/ao) are shown 
below^. The covariance C is chosen such that the samples u almost surely have 
regularity u G (D) for all s < where from left to right t = 2,1.5,1,0.5 

respectively. Here the samples are generated on [—1,1]^ 3 D and then restricted to 
D, for computational simplicity. 



3.1.2. Star-shaped prior. We now consider star-shaped inclusions, that is, inclusions 
parametrised by their centre and a radial function. These were studied in two- 
dimensions in the paper [8] to parametrise domains for a Bayesian inverse shape 
scattering problem. In [8] the authors prove well-posedness of the inverse problem 
in an infinite dimensional setting through the use of shape derivatives and Riesz- 
Fredholm theory. 

Let D C and Rd-i = (—7r,7r] x [0,7r]^“^ C Let h : ^ Rd-i be 

the continuous function representing the mapping from Gartesian to angular polar 
coordinates. Define the mapping A : Cp{Rd-i) x D ^ >S(D) by 

A{r,xo) = {x e D \ \x — xo\ < r{h{x — tq))} 

where Cp{Rd-i) is the space of continuous periodic functions on Rd-i- Then 
H(r, xo) describes the set of points in D which lie within the closed surface 
parametrised in polar coordinates centred at xq by 

r(0) = (0,r(©)), OeRd-i. 

^Given a measure fi on (AT, A’) and a measurable map F : (AT, A) —)► (Y, y) between measurable 
spaces, F^{/j,) denotes the pushforward of p, by F, i.e. the measure on (Y, y) given by F^{/j,){A) = 
/r(F-i(A)) for all A G A- 
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In two dimensions, we have Ri = (—tt, tt] and the mapping h : ^ i?i is given by 



where atan2 is the two-parameter inverse tangent function. 

In three dimensions, we have R 2 = (— tt, tt] x [0, tt] and the mapping h ^ R 2 
is given by 



Similar expressions for h exist in higher dimensions, though for applications we 
are only interested in the case d = 2, 3. 

Define now the map F 2 : Cp{Rd-i) x D ^ by 


F2{r,xo) = R+1a( r,xo) + U-^D\Ai r,xo) 
= {u+ - U_)lA(r,xo) +■“-> 


where > 0 are the scalar conductivity values. Again it can easily be seen 

that F 2 does indeed map into A{D). We claim that this map is continuous in the 
following sense: 

Proposition 3.2. Define the map F 2 : Cp{Rd-i) x D ^ as above. Let 

xq e D and let r G Cp{Rd-i) be Lipsehitz eontinuous. 

(i) Suppose that (re)e>o ^ Cp{Rd-i) is a sequenee of funetions sueh that 
Ike ~ ^lloo ^ 0- Then F 2 (r£,xo) ^ F 2 (r, tq) in measure^. 

(a) Suppose that (xo)e>o F D is a sequenee of points sueh that \xq — xo\ ^ 0. 
Then F 2 (r, Xq) ^ F 2 (r, tq) in measure. 

(Hi) Let {rs)e>o, (^o)e>o above. Then F 2 (re,Xo) ^ F 2 (r, xq) in measure. 

Proof. In order to show that a sequence of functions f^ : D ^ R converges to 
/ : D ^ M in measure, it suffices to show that there exists a sequence of sets 
C D with \Zs\ -A- 0 such that |/e — /| < Clz^- Then for each d > 0 we have 


\{x eD \ \f^{x) - f{x)\ > d}| <\{xe D \ \f^{x) - f{x)\ ^ 0}| < |Ze| ^ 0. 


(i) Fix the centre xq G D. Denote A{r) = A(r, xq). Let r G Cp{Rd-i) and let 
{'^s)s>o F Cp{Rd-i) be a sequence of functions such that Ik^ — r||oo ^ 0- 
Then there exists 7(5) ^0 such that |ke “'^||c» < 7(^)- By definition we then 
have 

r{x) — 7(5) < rs{x) < r{x) + 7(5) for all x G D and e > 0. 

It follows that we have the inclusions 


A{r - 7(5)) C A{r^) C A{r + 7(5)), 
A{r — 7(5)) C A{r) C A{r + ^{s)). 


Let A denote the symmetric difference. We deduce that 
A{rs)AA{r) C A(r + 7(^)) \ A{r — 7(5)). 


^A sequence of functions (fe)e>o, A : D ^ M, is said to converge in measure to a function 


/ : D —M if for all (5 > 0, \{x e D \ \fe{x) — f{x)\ > (5}| —0. Here \B\ denotes the Lebesgue 
measure of a set B C 
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Now the right-hand side is given by 
A{r + 'y{s))\A{r-j{e)) 

= {x e D \ r{h{x - xq)) - 7 (e) < \x - xo\ < r{h{x - xq)) + 7 (e)} • 

As 5 ^ 0, this set decreases to the boundary set 

dA{r) = {x e D \ \x — xo\ = r{h{x — tq))} . 

Since the graph of a continuous function has Lebesgue measure zero, we deduce 
that \dA{r)\ = 0. It follows that 

lim \ A{rs)AA{r) \ = 0. 

To conclude, note that 

|-^2(^e5 ^o) -^ 2(^5 ^ 0 ) I ^ 1"^+ I |I-A(r£) \ A{re)A>.A{r) • 

(ii) Let r G Cp{Rd-i) be Lipschitz continuous. Denote A{xq) = A(r, tq). Let 
(xq) C D be a sequence of points such that \xq — xq\ ^0. Note that we may 
write 

^(^ 0 ) = {x e D \ \x — xl\< r{h{x - xg))} 

= {x \ \x — Xq\ < r{h{x — Xq))} n D 
= {{xq - Xq) {x \ \x - Xo\ < r{h{x - To))}) n D 

=: {{xq - To) + A(to)*) n D. 

By the distributivity of intersection over symmetric difference, we then have 
that 

A{xq)AA{xo) = [((Tg - To) + A(to)*) n D]A[A{xoy n D] 

= [((Tg - To) + A(to)*)AA(to)*] n D 
C {{xq - To) + A(to)*)AA(to)*. 

Therefore, using Theorem 1 from [39], we see that 

|A(Tg)AA(To)| < |((Tg - To) + A(to)*)AA(to)*| 

< \xl - xo\'H'^~^{dA{xoy) 

where is the {d — l)-dimensional Hausdorff measure. Since we assume 

that r is Lipschitz, the surface area H^“^(9 A(to)*) of the boundary of A(to)* 
is finite, and so it follows that 

lirn |A(to)AA(to)| = 0. 

As before, we conclude by noting that 

|i^ 2 (r,Tg) -F 2 (r,To)| < |r+ -U-WlAixf^) -tAixo)\ =CtAixf^)AAixo)‘ 

(iii) We have that 

|i^2(r£,T§) -F2(r,To)| < |F2(r^,T§) -F2(r,T§)| + |F2(r,T§) -F2(r,To)| 

— ^{^A(r£,XQ)AA(r,XQ) T ^A{r^XQ)AA{r,xo)) 

— [A(rs ,XQ)AA(r,XQ)]U[A(r,XQ)AA(r,xo)] • 

Now note that 

|A(r£,^o)AA(r,^o)| < |A(r£,^o)*AA(r,^o)*|- 
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The right hand-side is independent of by translation invariance of the 
Lebesgue measure. By the same argument as part (i) we conclude that it 
tends to zero. We then have that 

\[A{re,xl)AA{r, xl)] U [A{r, a:g)AA(r, a;o)]| 

< \M're,xl)AA{r,xl) \ + \A{r,xl)AA{r,xo)\ 

< sup |^(re, 2 /o)A^(r, 2 /o)| + |-4(r,a;o)AA(r,a;o)| 

voeD 

which tends to zero by the discussion above and part (ii). 

□ 


Remark 3.3. Above we assumed that r : Rd-i M was Lipschitz continuous. 
This assumption is only used in the proof of part (ii) of the proposition. If the 
centre of the star-shaped region is known, this assumption may then be dropped to 
allow for rougher boundaries. 

We need to choose a prior measure yo on (r, xq). We equip Cp{R(i-i)xD with the 
norm ||(r, to)|| = max{||r||oo 5 |^o|} and corresponding Borel cr-algebra. We assume 
that r and Xq are independent under the prior so that we may factor /ig = ctq (8) tq 
where erg is a measure on Cp{Rd-i) and rg is a measure on D. We will assume that 
cTg is such that erg(5) > 0 for all balls B C Cp{Rd-i)‘ 

Example 3.4. Consider the case D = C Suppose that r is drawn from a 

log-Gaussian measure erg on Cp((— tt, tt]), and xq is drawn from rg = I/([—0.5, 0.5]^). 
Note that [—0.5, 0.5]^ C 5(0,1). Typical samples from (yo) are shown below. 
The covariance of <jg is chosen such that the samples r almost surely have regularity 
r G ((— TT, tt]) for all s < t, where from left to right t = 2.5,2,1.5,1 

respectively. 








w 




3.1.3. Level set prior. We finally consider the case where the inclusions can be 
described by a single level set function, as in [17]. Let n G N and fix constants 
—oc = cg < Cl < ... < = oc. Given r : 5 ^ M, define Di C D hy 

Di = {x £ D \ Ci-i < u{x) < Ci}, i = 1,..., n 

so that D = UlLi n Dj = 0 for i ^ j > 1 . Define also the level sets 

= Di0 Dipl = {x G 5 I u{x) = Ci}, i = 1,..., n — 1. 

Now given strictly positive functions /i,..., /n G C^{D), we define the map F 3 : 
C^(D) ^ A{D) by 

n 
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Since each / is continuous and strictly positive on a compact set they are uni¬ 
formly bounded above and below by positive constants, and so F 3 does indeed map 
into A{D). 

In this paper we are primarily concerned with the case of binary fields, n = 2 
and fi constant above, however the theory in proved in the general case. We have 
the following result regarding continuity of this map, by the same arguments as in 

[ 17 ]. 

Proposition 3.5. Define the map F 3 : C^{D) A{D) as above. Let u G C^{D) 
he sueh that \D^\ = 0 for i = 1 ,... ,n — 1 . Suppose that (us)s>o ^ C^{D) is an 
approximating sequenee of funetions so that — r||oo ^0. Then F^{uA) F^{u) 
in measure. 

Proof. Denote by and ^ the sets as defined above associated with the ap¬ 
proximating functions u^. We can write 

n n n 

Fs{Ue) - F 3 {u) = 'T.'P.ifi - = 'P. ifi - fj)'^Di,^nDj ■ 

i=l j=l hj=l 

Since \\us — r||oo ^ 0, there exists 7 ( 5 ) ^ 0 with \\us — r||oo < 7 (^)- Then we have 
for all X G D and 5 > 0 

u{x) — 7 ( 5 ) < Us{x) < u(x) + 7 (^). 

Hence for \j — i\> 1 and e sufficiently small, H Dj = 0. If jj — = 1 , then 

A,£n A+i c A,£ ■■= {x g D\ci< u{x) <Ci + 7(e)} D°, 

Di^s n A-i C A-i,£ ■■= {x G D\ci- 7(e) < u{x) < c*} 0. 

By the uniform boundedness of the (/i), for sufficiently small 5 we can then write 

n— 1 n 

\Fs{uV - F 3 {u)\ /i+il^Ae + E I/* - 

1=1 ’ i =2 

( 6 ) < Clz^ 

where 0 D is given by 

/n— 1 \ ^ \ 71—1 

A = U A,, u MJ A-i,e ^ U A°. 

\i=l / \i=2 / i=l 

By the assumption that \D^\ = 0 for all i, it follows that \Z^\ ^0, and so the result 
follows from the comment at the start of the proof of Proposition 3.2. □ 

Note that bound ( 6 ) actually above implies the slightly stronger result that, 
when the c^-level sets of u e X have zero measure, F 3 is continuous into Lp{D), 
1 < p < 00 , at A. The assumption that the level sets have zero measure is an 
important one, as illustrated by Figure 2 : an arbitrarily small perturbation of u 
can lead to an order 1 change in Fs{u). 

In the Bayesian approach we are taking to this problem, we may choose a prior 
measure on u such that, almost surely, the Lebesgue measure of the level sets is 
zero. This is shown to hold for non-degenerate Gaussian measures in [17]. As a 
result, F 3 will be almost surely continuous under the prior, and this is enough to 
give the measurability required in Bayes’ theorem, as shown in [17]. 
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Figure 2. The discontinuity of F3 into Lp{D). 


As in the log-Gaussian case, we take our prior measure /io on r to be a Gaussian 
measure N{mQ^Co) 011 C^{D). 

Example 3.6. Gonsider the case D = 5(0,1) C n = 2, ci = 0, /i = 1 
and /2 = 2. Suppose that u is drawn from a centred Gaussian measure /io = 
A^(0,C) on C^{D). The covariance C is chosen such that the samples u almost 
surely have regularity u G (5) for all s < where from left to right t = 

4, 3, 2,1 respectively. As in the log-Gaussian case, here the samples are generated 
on [—1,1]^ 3 D and then restricted to 5, for computational simplicity. 



Remark 3.7. In the cases of the star-shaped and level set priors, we have assumed 
the values that the conductivity takes are known a priori. This may be appropri¬ 
ate in certain situations, for example in medical contexts where conductivities of 
different types of tissue are known [7]. However it may be preferable to have the 
flexibility of treating the conductivity values as part of the inference. The theory 
does not become significantly more involved, as it can be shown that the conclu¬ 
sions of propositions 3.2 and 3.5 still hold when we allow for perturbation of the 
conductivity values as well. We work with the case of fixed conductivity values for 
clarity of presentation, but we note that it is possible to place a prior upon these. 

3.2. The Likelihood and Posterior Distribution. The inverse problem was 
introduced at the beginning of the section. Now that we have introduced prior 
distributions, we may provide the Bayesian formulation of the problem. 

Let X be a separable Banach space and F : X —> A{D) a map from the space 
X where the unknown parameters live to the conductivity space. Choose a set of 
current stimulation patterns G j = 1,..., J and let Mj : A{D) H denote 
the solution map when using stimulation pattern Recall the projection map 

n : H ^ was defined by n('i;, V) = V. 
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The data yj from the jth stimulation pattern is assumed to arise from the map 
Qj : X ^ Qj = Ho Xij o F, via 

Vj = Sj{u) ^r]j, r]j - A^( 0 ,ro) iid. 

We concatenate these observations to get data y G given by 

y = Q{u) + 7^, y ^Qo N{0, T) 

where T = diag(ro,..., Tq) and Q \ X ^ This coincides with the setup at the 

start of the section, with a = F{u). 

Assume that ix ^ /io, where /xq is independent of Qo- From the above, we see 
that ^|xx - Qn := N{Q{u)^ r). We use this to find the distribution of xx|^. First note 
that 

S(^)=exp(-^(«; 2 /) + ^Hr) 

where the potential (or negative log-likelihood) ^ : X x T ^ M is given by 

( 7 ) ^{u-,y) = ^\Q{u)-y\l := {g{u) - y)\'^. 

Then under suitable regularity conditions, Bayes’ theorem tells us that the distri¬ 
bution /x^ of XXI ^ is as given below: 

Theorem 3.8 (Existence and Well-Posedness). Let (X, F, /xq) denote any of the 
probability spaees assoeiated with any of the three priors introdueed in the previous 
subseetion, and let ^ : X x Y R be the potential (7) assoeiated with the eorre- 
sponding forward map. Then the posterior distribution pT of the state u given data 
y is well-defined. Furthermore, jF <C /xq with Radon-Nikodym derivative 

(8) exp(-$(M; y)) 

where for y Q^^-a.s., 

Zfx := / exp(-$(M; 2 /))Mo(dM) > 0. 

Jx 

Additionally, the posterior measure jF is loeally Lipsehitz with respeet to y, in the 
Hellinger distanee: for all y,y' G Y with max{|^|r, |^^|r} < P, there exists C = 
C{p) >0 sueh that 

d}ie\\{p^F^ ) < Ct\y - y'\r. 


In the proof of the above theorem we will make use of the following version of 
Bayes’ theorem from [10]. 


Proposition 3.9 (Bayes’ theorem). Define the measure zxo(dxx, dy) = po{du)Qo{dy) 
on X X Y. Assume that ^ : X x Y ^ M UQ-measurable and that, for y Q-a.5. 

Z^= exp(-<h(xx;^))/xo(dxx) > 0. 

Jx 

Then the eonditional distribution of u\y exists and is denoted by pA. Furthermore 
pA po and, for y Qo-a.s., 


dpy 

dpo 


1 


exp(-$(xx;^)) . 


We need to verify that the assumptions of this theorem are satisfied. To proceed 
we first give some regularity properties of the potential 
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Proposition 3.10. Let (X, /io) denote any of the probability spaees assoeiated 
with the priors introdueed in the previous subseetion. Then the potential ^ \ XxY 
M assoeiated with the eorresponding forward map, given by (7), admits the following 
properties. 

(i) There is a eontinuous K : IR+ x IR+ ^ IR+ sueh that for every p > 0, u ^ X 
and y with \y\r < p, 

0 < ^{u\y) < K{p, ||m||x)- 

In the eases F = F 2 and F = Fs, K has no dependenee on ||r||x- 
(a) For any fixed y ^Y, ^{';y) : X ^ M eontinuous po-almost surely on the 
probability spaee {X,F,po). 

(Hi) There exists C : x M+ ^ M+ sueh that for every yi,y 2 ^ T with 

max{|^i|r, |^ 2 |r} < P, ctnd every u e X, 

\^{u;yi) - $(«; 2 / 2)1 < C{p, ||■u|U)|2/l - 2/2|r- 

Moreover, C{p, || • ||x) e L^^iX) for all p>0. 

Proof (i) From equation (5) in the proof of Proposition 2.7, we see that there 
exists C > 0 such that each Mj : A{D) H satisfies 

\\Mj{ai) - Mj{a2)\U < C||Mj((T2)||*lki - cr2||oo 

for all cTi, (72 G A{D). Taking (72 = 1, say, we deduce that 

||A^i(cTi)|U < C\\M,il)\U\\a, - llloo + ||A^i{l)|U < C(1 + ||ai|U). 

Hence ||(7||oo < P implies that ||Xlj((7)||* < C{1 F p). By Corollary 2.8, it 
follows that n o Xij : A{D) is bounded on bounded sets with respect 

to II • Iloo for all j. 

In the case F = Fi, if r G X then ||F(r)||oo < It follows that 

|0(R)|r < maxj \Qj{u)\T < C{1 + 

Now note that 

^{u‘y)<\g{u)\lF\y\l 
Then for any y ^Y with |^| < p, we may bound 

i>(w; y) < 0(1 + + P) =: K{p, ||ti||x). 

In the cases F = F 2 and F = F 3 , we have that ||F(r)||oo is bounded uniformly 
over R G X and so \g{u)\r < maxj |0j('u)|r < C. Hence we obtain the bound 

<^>{u-,y)<C{l+p^)=:K{p). 

(ii) Let u po and suppose F : X ^ A{D) is such that \\us — u\\x —> 0 implies 
that F{uA) ^('^) either uniformly or in measure. Then Proposition 2.7 
tells us that Mj o F : X ^ H is continuous at u for each j. The projection 
H : H ^ is continuous, and so Qj = H o Mj o F is continuous at u for 
each j. In §3.1 it is shown that this is true for F = Fi and F = F 2 for any 
u. For F = F 3 it is only true at points u whose level sets have zero measure, 
however since we are assuming r ^ po, a Gaussian measure, it follows from 
Proposition 7.2 in [17] that u po-almost surely has this property. 
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(iii) Let u e X and ^1,^2 ^ ^ with max{|^i|r, |^2|r} < P- Then we have 

\^{u-,yx) - ^{u\y2)\ = ^1(2/1+2/2- 2 ^('u),2/1-2/2)r| 

^ 2^\yi\T + \y 2 \T + 2|^(M)|r)|2/i - 2/2|r 
< (p + \Q{u)\T)\yi - 2 / 2 |r 
=: C{p,u)\yi - 2/2|r 

We now consider cases separately based on the prior. In the log-Gaussian case, 
we may bound 


«) < ca + p + e'HU) =: C(p, IMIx) 

using the bound from the proof of part (i). Square-integrability of C(p, || • ||x) 
follows since Gaussians have exponential moments. 

In the star-shaped and level set prior cases, we have that \Q{u) \ is bounded 
uniformly by a constant. We may hence bound C{p,u) above by some 
C{p,\\u\\x) •= C(1 + p) that is independent of ix, and so again the square- 
integrability follows. 

□ 


Proof of Theorem 3.8. Define the product measure = yo{du)Qo{dy) on 

X xY. We showed in Proposition 3.10 that ^{';y) : X ^ M is almost-surely 
continuous under the prior for all y e and •) : Y ^ M is locally Lipschitz 
for all XX G X. Together these imply that ^ : X x T ^ M is almost-surely jointly 
continuous under uq. To see this, let {u,y) £ X x Y and let (xXn,Pn)n>i ^ X x T 
be an approximating sequence so that ||xxn — u\\ X ^0 and \yn — y\r 0. Then we 
have 

\^{Un,yn) - ^(W,2/)| < \^{Un,yn) “ ^{Un,y)\ + \^{Un,y) “ 

The second term tends to zero /xo-almost surely by continuity. For the first term, 
note that the sequences (||xXn||A)n>i and (|pn|r)n>i are bounded, by K and R 
respectively, say. Then we can use the local Lipschitz property to deduce that 

\^{un,yn) - ^{Un,y)\ < C{R,K)\yn - y\r 

since (7(', •) : MxM ^ M, as defined in the proof of Proposition 3.10, is monotonically 
increasing in both components. Therefore this term tends to zero, and we obtain 
the desired continuity. It follows, see for example Lemma 6.1 in [17], that $ is 
z/q- measurable. 

For the lower bound on we consider cases separately based on the prior. First 
we consider the log-Gaussian and level set prior cases so that /xq is Gaussian. Let 
B C X he any ball. Fix any p > \y\r and define 

R = sup if (/o, ||u||x) 

ueB 
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where K is the upper bound from Proposition 3.10(i). This supremum is finite by 
the continuity of K. Then we have 

/ exp(-$(u;y))/io(dw) > / exp(-$(u; y))/io(du) 

Jx J B 

> f exp{-K{p,\\u\\)) noidu) 

Jb 

> exTp{-R)po{B). 

Since po is Gaussian, po{B) > 0 and so > 0. 

In the star-shaped prior case, proceed as above but take B = Bi x D where 
Bi C C°p{Rd-i) is any ball. Then we have 

yo(-B) = (o-Q X to){Bi X D) = cfq{Bi)tq{D) > 0 

by the assumption that do assigns positive mass to balls, and so again > 0. The 
above hold for all ^ G T, and so in particular for y Qo-almost-surely. We may now 
apply Bayes’ Theorem 3.9 to obtain the existence of . 

The proof of well-posedness is almost identical to that of the analogous result 
Theorem 2.2 in [17] and is hence omitted. □ 

4. Numerical Experiments 

We investigate the effect of the choice of prior on the recovery of certain binary 
conductivity fields. The specific fields we consider are shown in Figure 3, where 
blue represents a conductivity of 1 and yellow a conductivity of 2. Simulations are 
performed using the EIDORS software [1] to solve the forward model using a first 
order finite element method; a mesh of 43264 elements is used to create the data 
and a mesh of 10816 elements is used for simulations in order to avoid an inverse 
crime [20]. 

In subsection 4.1 we describe the MCMC sampling algorithm that we will use. 
In subsection 4.2 we define the parameters we will use for the forward model and 
the MCMC simulations. We also describe how the data is created, and define our 
choices of prior distributions. Finally in subsection 4.3 we present the results of the 
simulation, looking at quality of reconstruction, convergence of the algorithm and 
some properties of the posterior distribution. 

4.1. Sampling Algorithm. We aim to produce a sequence of samples from on 
X, where is given by (8). We make use of the preconditioned Crank-Nicolson 
Markov Chain Monte Carlo (pCN-MCMC) method. The pCN-MCMC method is 
a modification of the standard Random Walk Metropolis (RWM) MCMC method 
which is well-adapted to Gaussian priors in high dimensions. It was introduced in [5] , 
and its dimension independent properties are analysed and illustrated numerically in 
[14] and [9] respectively; the pCN nomenclature was introduced in [9]. In the case of 
the star-shaped prior, we use a Metropolis-within-Gibbs algorithm [43], alternately 
updating the field with the pGN method above and updating the centre with the 
standard RWM method. 

An advantage of these MGMG methods is that derivatives of the forward map 
are not needed, only black-box solution of the forward model. However in order 
to accurately compute some quantity of interest, such as the conditional mean, we 
may need to produce a very large number of samples and tuning the algorithm to 
minimise this effect is important. For this reason we compute the effective sample 
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(a) Conductivity A (b) Conductivity B 

Figure 3. The two binary fields we attempt to recover. Con¬ 
ductivity A is drawn from the star-shaped prior, with ao = 
[A^(0.5, 10^ • (30^ — Ad)~^)] , h{z) = (1 + tanhz)/ 2 , and tq = 
[/([—0.5, 0.5]^). Conductivity B is constructed explicitly, rather 
than being drawn from a prior. 


size from the integrated autocorrelation (neglecting a burn-in period) of a quantity 
of interest, as in [ 21 ]. 


4.2. Data and Parameters. We work on a circular domain of radius 1, with 16 
equally spaced electrodes on its boundary providing 50% coverage. We take all 
contact impedances zi = 0.01. We stimulate adjacent electrodes with a current of 
0.1, so that the matrix of stimulation patterns I = G is given by 


/ +1 

-1 


/ = 0.1 X 


0 

+1 


0 -1 


V 0 


0 0 


0 \ 

0 

0 

+1 

- 1 J 


The conductivity is chosen such that it takes values 1 and 2 . We perturb the 
measurements with white noise r] ^ A^( 0 , 7 ^/), 7 = 0 . 0002 , so that the mean relative 
error on both sets of data is approximately 10%. The true conductivity fields used 
to generate the data, henceforth referred to as Conductivity A and Conductivity 
are shown in Figure 3. In all cases we generate N = 2.5 x 10^ samples with a 
burn-in of /cq = 5 x 10 ^ samples. 

Our priors on fields will make use of Gaussians with covariances of the form 
(9) C = q{r^-A)-^. 


These are essentially rescaled Whittle-Matern covariances [31], with r representing 
the inverse length scale of the samples, a proportional to their regularity, and q 
proportional to their amplitude. 

In what follows, denote by the Laplacian with Neumann boundary conditions 
on [— 1 , 1 ]^, restricted to D, so that its domain is given by 


B{An) = \ u\d 


uGHy[-i,in^ = o 
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Defining the Laplacian first on a square and then restricting to D will allow for 
efficient generation of Gaussian samples via the fast Fourier transform. Note that 
if we were to consider priors of the form (9) with r = 0, we should restrict V{An) 
further to ensure the invertibility of An- 

Additionally, denote hy Ad the Laplacian with Dirichlet boundary conditions on 
Ri = (—7r,7r], so that its domain is given by 

= {r G iL^((—7r,7r]) | u{—7t) = u{7r) = O}. 

4.2.1. Gaussian prior. States are defined on a grid of 2^ x 2^ points. For both 
simulations the pCN jump parameter P is taken to be 0.01, with choice of prior 

/io = exp^ [A^(0.51og2,10^^ • (40^ - • 

4.2.2. Star-shaped prior. Radial states are defined on a grid of 2^ points. For Con¬ 
ductivity A, we choose the pCN jump parameter (3 = 0.03 and the RWM jump 
parameter 5 = 0.01. For Conductivity B we choose p = 0.01 and 5 = 0.005. For 
both simulations we use the choice of prior /ig = ctq x tq, with 

(TO = h* [A^(0.5,10^ • (30^ - Ad)-'')] , tq = C/([-0.5,0 .S]^), 
where h{z) = (1 +tanh2:)/2. 

4.2.3. Level set prior. States are defined on a grid of 2^ x 2^ points. For both 
simulations the pCN jump parameter P is taken to be 0.005, with choice of prior 

^^o = N{0, {35 ^-An)-^). 

4.3. Results. 

4.3.1. Reeovery. Figures 4 and 5 show conductivities arising from the MCMC chains, 
and Figure 6 shows the values of the misfit at the different sample means. We 
consider two different types of mean conductivities. First the sample means are 
calculated in the sample spaces Xi and then pushed forward to the conductivity 
space by the maps Fi^ so that we produce estimates of F^(E(r)). This preserves 
the binary nature of the fields in the cases of the star-shaped and level set priors. 
We also consider the sample means of the pushforwards of the posteriors by 
E{Fi{u)), which do not preserve the binary property but provide some visualisation 
of the interface uncertainty in the posterior. 

For Conductivity A, the sample mean F2 (E(r)) arising from the star-shaped 
prior provides a better reconstruction than the other two prior choices. This is 
expected, since the true conductivity was drawn from this prior. Whilst the sample 
means arising from the level set prior are fairly close to the true conductivity (both 
visually and in terms of 4>), the boundary of the interface appears to have too large 
a length-scale. Appropriate choice of prior length-scale is a key issue the with the 
level set method; treating the length-scale hierarchically as another unknown in the 
problem may be beneficial. On the other hand 

The sample means arising from the Gaussian prior fail to recover both the sharp 
interface and the values of the conductivity, which is reflected in the large values 
4> takes. Nonetheless, general qualitative properties of the true conductivity can 
be seen in the sample means. As with the level set prior, reconstruction could be 
improved by treating hierarchically the parameters in the prior such as length-scale, 
marginal variance and mean. 

For Conductivity B, the level set prior is most effective in the reconstruction, 
since a specific number of inclusions isn’t fixed a priori as it is for the star-shaped 
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prior. Again the Gaussian prior fails to recover both the sharp interface and the 
values of the conductivity, however it appears to do a better job than the star-shaped 
prior at identifying the location and shape of the two inclusions. 

In both of the above cases, even though the individual samples coming from 
using the level set prior contain many small inclusions, these do not show up in the 
sample means F3 (E(r)), though they are visible in the means E{Fs{u)). 

In terms of the values of the misfit as shown in Figure 6, in all but one of the 
above cases the mean E(Fi(R)) provides a better fit than Fi{E{u)) despite the lack 
of the binary property. The quality of the recovery cannot be assessed solely on 
the misfit however; an advantage of the Bayesian approach is that both means are 
available (along with the full posterior), and whichever is more appropriate could 
depend on the context. 

4.3.2. Convergence. In Figure 7, we show the approximate effective sample size 
(ESS) associated with different quantities of interest. For all choices of prior, these 
are significantly smaller than the total 2.5 x 10^ samples generated. Many more 
samples may hence be required to produce accurate approximations of the posterior 
mean. 

The chain associated with the star-shaped prior results in the largest ESS, likely 
because we are only attempting to infer 2^+2 parameters rather than 2^^ parameters 
as in the log-Gaussian and level set cases. 

In order to accelerate the convergence of the MGMG we can adjust the jump 
parameters (3 and 5. Larger choices of these parameters mean that accepted states 
will be less correlated with the current state, however the proposed states are less 
likely to be accepted. The choice /3 = 1 in pCN produces proposed states that are 
independent of the current state, but dependent upon how far the prior is from the 
posterior, very few or no states may be accepted so that the chain never moves. 
Similarly, smaller choices of these jump parameters mean that more proposals will 
be accepted, but the states will be more correlated. A balance hence must be 
achieved - in our simulations we choose the parameters such that approximately 
20-30% of proposals are accepted, though in general the optimal acceptance rate is 
not known [3]. 

Alternatively, reconstruction may be accelerated by looking at an approximation 
of the posterior instead of the exact posterior, for example using the ensemble 
Kalman filter [16] or a sequential Monte Carlo method [4]. We could also initialise 
the MGMG chains from EnKE estimates to significantly reduce the burn-in period 
and hence computational cost. If the derivative of the forward map is available. 
Hybrid Monte Carlo (HMC) methods could be used to accelerate the convergence 
[11]. Emulators could also be used to reduce the computational burden of derivative 
calculation, allowing the use of geometric MGMG methods such as Riemannian 
Manifold Hamiltonian Monte Carlo (RHMC) and Lagrangian Monte Carlo (LMC) 
[27]. 

4.3.3. Posterior Behaviour. In Eigures 8-10 we show kernel density estimates for a 
number of quantities associated with each posterior distribution. These are calcu¬ 
lated using subsequences of the MGMG chains, with lengths of the same order of 
magnitude as the effective sample sizes, in order to avoid over-fitting. The most 
regular densities arise in the star-shaped case, with the distribution of all quan¬ 
tities appearing to be very close to uni-modal. More irregularity is seen for the 
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log-Gaussian case, especially in the joint distributions, but they are still very close 
to uni-mo dal. 

The least regular, more multi-modal densities come from the level set prior. One 
reason for this is likely the lack of identifiability of the level set function: the forward 
model only ‘sees’ the zero level set of the state, and hence cannot distinguish between 
infinitely many different states. The prior can however distinguish between these 
states, and will weight them appropriately, which can help explain the shape of the 
posterior densities. 


5. Conclusions 

The primary contributions of this paper are: 

• We have formulated the EIT problem rigorously in the infinite dimensional 
Bayesian framework. 

• We have studied three different prior models, each with their own advan¬ 
tages and disadvantages based on prior knowledge and the nature of the 
field we are trying to recover. 

• With each of these choices of prior we obtain well-posedness of the problem. 
We can obtain well-posedness using additional prior models, as long as the 
mapping from the state space to the conductivity space has appropriate 
regularity. 

• The infinite dimensional formulation of the problem leads to the use of 
state of the art function space MCMC methods for sampling the posterior 
distribution. 

• Simulations performed using these methods illustrate that the conditional 
mean provides a reasonable reconstruction of the conductivity, even with 
fairly significant noise on the measurements. They also illustrate the fact 
that the choice of prior has a significant impact on reconstruction and, 
in particular, that the geometric priors (star-shaped and level set) can be 
particularly effective for the (approximately) piecewise constant fields that 
arise in many applications. 

Future research directions could include the following: 

• Sampling the exact posterior distribution using MCMC can be computa¬ 
tionally expensive. Methods that approximate the posterior may be as 
effective for calculating quantities such as the conditional mean, with much 
lower computational load. The relative effectiveness versus cost of different 
methods could be studied. This could be especially important for simula¬ 
tions in three dimensions, where forward model evaluations are even more 
expensive. 

• When using the level set prior, the length scale of samples could be treated 
hierarchically as an additional unknown in the problem. 

• The star-shaped prior could be extended to describe multiple inclusions, 
either with the number of inclusions fixed or as an additional unknown. 

• The Calderon problem could be considered in place of the CEM. The 
data could then be either the full Neumann-to-Dirichlet (or Dirichlet-to- 
Neumann) map, or a finite number of pairs of Dirichlet and Neumann 
boundary values. In order to perturb the measurements with noise, the 
former case would require the definition of Gaussian distributions on spaces 
of operators, and the latter Gaussian distributions on manifolds. 
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(a) (Top) Estimates of Fi(E(i^)). (Bottom) Estimates of E(Fi(i^)). 



(b) Posterior samples. 

Figure 4. Recovery of Conductivity A. From left to right, the 
log-Gaussian, star-shaped and level set priors are used. 
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(a) (Top) Estimates of Fi(E(i^)). (Bottom) Estimates of E(Fi(i^)). 



(b) Posterior samples. 

Figure 5. Recovery of Conductivity B. From left to right, the 
log-Gaussian, star-shaped and level set priors are used. 
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Conductivity A 

Conductivity B 

Log-Gaussian prior 

115960 

122180 

Star-shaped prior 

228.84 

45013 

Level set prior 

284.28 

193.03 


(a) The values of $ at Fi(E(u)). 



Conductivity A 

Conductivity B 

Log-Gaussian prior 

111960 

117950 

Star-shaped prior 

284.44 

41408 

Level set prior 

137.93 

124.42 


(b) The values of $ at E(Fi(u))). 


Figure 6. The values of the misfit at the different mean con¬ 
ductivities, for the different conductivities and prior distributions. 


Quantity 

Estimated ESS 

Quantity 

Estimated ESS 

u{0,l) 

40.0 


241.7 

u{0,2) 

90.7 

^(2) 

89.6 

w(l,l) 

35.4 

r(l) 

101.1 

u{l,2) 

44.5 

r(2) 

179.4 

u{l,3) 

36.0 

r(3) 

277.8 

u{2,l) 

101.9 

f(4) 

214.8 

u{2,2) 

37.9 

f(5) 

146.7 

u{2,3) 

89.7 

r(6) 

146.7 


(a) Log-Gaussian prior (b) Star-shaped prior 


Quantity 

Estimated ESS 

m(0 , 1) 

26.4 

u{0,2) 

28.9 

6 (1,1) 

27.2 

6 (1,2) 

23.5 

6(1,3) 

23.5 

6 (2,1) 

26.0 

6 (2,2) 

26.6 

6(2,3) 

24.1 


(c) Level set prior 


Figure 7. (Conductivity B) The estimated effective sample size 
for each chain, approximated using a variety of quantities, for the 
different choices of prior. In all cases 2.5 x 10^ total MCMC samples 
are produced, with the initial 5 x 10^ discarded as burn-in. 
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Figure 8. (Conductivity B, log-Gaussian prior) Kernel density 
estimates associated with six Fourier coefficients of u. The diag¬ 
onal displays the marginal densities of each coefficient, and the 
off-diagonals the marginal densities of corresponding pairs of coef¬ 
ficients. 
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Figure 9. (Conductivity B, star-shaped prior) Kernel density es¬ 
timates associated with the centre and four Fourier co¬ 

efficients of r. The diagonal displays the marginal densities of 
each quantity, and the off-diagonals the marginal densities of cor¬ 
responding pairs of quantities. 
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Figure 10. (Conductivity B, level set prior) Kernel density es¬ 
timates associated with six Fourier coefficients of u. The diago¬ 
nal displays the marginal density of each coefficient, and the off- 
diagonals the marginal densities of corresponding pairs of coeffi¬ 
cients. Axes are rescaled by 10^ for clarity. 
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